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Phase diagrams of binary mixtures of oppositely charged colloids are calculated theoretically. 
The proposed mean-field-like formalism interpolates between the limits of a hard-sphere system at 
high temperatures and the colloidal crystals which minimize Madelung-like energy sums at low tem- 
peratures. Comparison with computer simulations of an equimolar mixture of oppositely charged, 
equally sized spheres indicate semi-quantitative accuracy of the proposed formalism. We calculate 
global phase diagrams of binary mixtures of equally sized spheres with opposite charges and equal 
charge magnitude in terms of temperature, pressure, and composition. The influence of the screen- 
ing of the Coulomb interaction upon the topology of the phase diagram is discussed. Insight into 
the topology of the global phase diagram as a function of the system parameters leads to predictions 
on the preparation conditions for specific binary colloidal crystals. 

PACS numbers: 64.70. pv, 64.70.D-, 64.70.K- 



I. INTRODUCTION 



Colloidal suspensions are interesting experimental re- 
alizations of many-particle systems because, in contrast 
to molecules, they are directly accessible to observations 
and the interactions can be tuned in a wide range by 
means of the preparation procedure and the experimen- 
tal conditions |lH4|. The flexibility of colloidal suspen- 
sions, however, implies many system parameters to be 
specified, such as size and shape of the colloids as well 
as strength and functional form of the interaction po- 
tential [5|-|8( . In addition, the number of parameters in- 
creases enormously if one considers mixtures of different 
colloidal species instead of a pure system of only one type 
of colloid. Moreover, the dimension of the phase diagram 
increases in parallel with the number of colloidal species 
in the suspension. The high dimensionality of the pa- 
rameter space and the phase diagrams call for a simple 
theory to get an overview of the global phase behavior of 
colloidal suspensions because a systematic scan by means 
of experiments or computer simulations is precluded. 

In the present work such a simple theory is formulated 
for the case of binary mixtures of spherical colloids of the 
same radius and opposite charges of equal magnitude. 
The global phase diagrams for this setting are in terms 



of the temperature, the pressure, and the composition 
and there is only one additional parameter describing the 
screening of the Coulomb interaction. The reason to con- 
sider this comparatively low-dimensional problem here is 
that the reliability of the approach should be assessed 
with respect to available results on special cases obtained 
by means of computer simulations [9|, [l0( and Madelung 
sum calculations p], LLD, |T2j ■ An extension of the method 
to multi-component, size- and charge-poly disperse col- 
loidal suspensions is readily achieved along the same 
lines. 

The proposed mean-field-like formalism correctly in- 
terpolates between the hard-sphere limit at high tem- 
peratures and a Madelung-type description at low tem- 
peratures. As a mean-field theory, the present approach 
underestimates fluctuations such that critical phenom- 
ena and quantitative aspects of phase transitions are cov- 
ered in lowest order only. The over-all topology of the 
phase diagram, in particular the types of phases present, 
however, is not expected to be influenced qualitatively. 
The simple and transparent formalism turns out to be 
in semi-quantitative agreement with the computer sim- 
ulation results of Ref. [l(|; this gives confidence that a 
possible extension to multi-component, size- and charge- 
polydisperse colloidal suspensions leads to a reliable de- 
scription. 

The model, the formalism, and further technical details 
are introduced in Sec.fll] Results are discussed in Sec. lIIII 
and Sec. II VI concludes with a summary. 
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II. FORMALISM 



Fluid phase 



A. Model 

Consider a binary mixture of N\ positively charged col- 
loids and N2 negatively charged colloids in a solvent of 
volume V, temperature T, Debye screening length k -1 , 
and dielectric constant e. For simplicity we restrict at- 
tention to equally sized spherical colloids of radius a with 
exactly opposite valency Z\ — Z > and Z2 = — Z < 
for the positively and negatively charged colloids, re- 
spectively. For later use we define the sign of the col- 
loidal charge <?,; = Zi/Z,i E {1,2}. The thermodynamic 
state can be characterized by the total volume fraction 
<j) = 4ira 3 N/(3V) and the mole fraction of the positively 
charged colloids x% = x = N\/N where N — Ni + N2 
is the total number of colloids; the mole fraction of the 
negatively charged colloids is given by xi — 1 — x. The 
osmotic pressure p will be reported in terms of the di- 
mensionless combination p* = 47rpa 3 /(3fcsT) where ks 
is the Boltzmann constant. The effective colloidal inter- 
actions are assumed to be pairwise additive with hard- 
core and screened- Coulomb contributions but without 
Van der Waals attractions; these conditions are realized 
in index-matched solvents. In terms of the dimensionless 
temperature T* = 4ire ekBTa(l + na) 2 /(Ze) 2 with the 
permittivity of the vacuum £0 and the elementary charge 
e, the pair potential can thus be written as 



PUij(r) = < q%qj exp(-/ta(r/q- 2)) 



, r < 2a 
,r>2a ' 



(1) 



T* r / a 

where /3 = l/(fcsT) is the inverse temperature. 

B. Gibbs free energy 

The phase behavior of the system under consideration 
will be calculated from the Gibbs free energy per particle 
per fcflT 



g(T*,p*,x) = f(T*,$(T*,p*,x),x) + 



V 



$(T*,p*,x) 



, (2) 



where f(T*, <f>, x) is the Helmholtz free energy per parti- 
cle per IzbT and $(T*,p*, x) is the total volume fraction 
for given (T*,p*,x), which is implicitly defined by means 
of the equation of state 



V 



•9f 



(3) 



ip=<S>(T* ,p* ,x) 



Below f(T*,(f>,x), and hence g(T*,p*,x) using Eqs. fl2J 
and ([3]), is calculated for (i) fluid, (ii) crystalline, and (hi) 
substitutionally disordered solid phases. The expressions 
involve a priori unknown constants which are fixed by 
fitting to the well-known hard-sphere fluid-crystal tran- 
sition. 



The Helmholtz free energy per particle of the fluid 
phase is approximated by the analytical solution of the 
mean spherical approximation (MSA) of hard-sphere 
Yukawa mixtures 11311 



f(T*,<f>,x) 



Cf + J^ Xi ln(xj) + ln(</>) 

i 

b(t) r 2 (r + 3Ka) 
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(4) 



with the integration constant c/, which will be fixed 
below, the excess internal energy due to the screened- 
Coulomb interaction B(T) given by Eq. (14) of Ref. |13| . 
and the effective screening strength T being the non- 
negative solution of the 6th-order algebraic equation (10) 
of Ref. [13j, which has to be solved numerically. Upon 
T* — >• 00 one finds r — ► such that Ecl (|4|) leads to the 
Carnahan-Starling equation of state [14j . 



2. Crystalline phases 

The Helmholtz free energy per particle of crystalline 
phases is approximated by the following mean-field func- 
tional of the one-particle density profiles Qi,i e {1,2}, 
where ft(r) is the number density of positively (i = 1) or 
negatively (i — 2) charged colloids at position r in space: 



/(T*,0)-/ hs (<^) 



(5) 



l - J&\ J dV ]T Qi&QiVWuW rj 

V V(r) z ' j 



= 2N ' d ' r 

V 



Here / hs is the hard-sphere Helmholtz free energy per 
particle per ksT in the crystalline phase, V is the system 
volume, and V(r) = {r' e V : |r' - r| > 2a}. Note that 
the concentration x of a crystal is fixed by the crystal 
structure and not a free parameter. Given the crystal 
structure C = C\ U C2 composed of the sublattices Ci,i € 
{1, 2} of positively (i — 1) and negatively (i = 2) charged 
colloids, an approximation to the density profiles is 



ft( r ) = J2$(r-s) 



(6) 



sect 



where S is the Dirac delta. The energetic contribution 
to the Helmholtz free energy per particle, given by the 
term on the right-hand side of Eq. ([5]), follows directly 
from Eqs. ([1]) and ((6|); it is the analog for the screened 
Coulo mb p otential Eq. (JT]) of the well-known Madelung 
energy [15j. The hard-sphere Helmholtz free energy per 
particle, / hs , is approximated by integrating the free- 
volume equation of state corresponding to crystal struc- 
ture C of closed-packed packing fraction <p cp [16[ with 
respect to the packing fraction </>. 



Finally the Helmholtz free energy per particle of crystal 
structure C reads 



f(T*,4>) 

c s + ln(» - 3 In ( 1 - I — 



(7) 



1/3' 



1 y^ V^ exp(— Ka(|s — s'l/a — 2)) 

+ W7F*^ qi<lj 4^ |s-s'|/a 

with the number of colloids per unit cell K c and the ab- 
breviation C* = {(s,s') € C t x Cj : s' € W c ,s ^ s'}, 
where Wc denotes the unit cell. The integration con- 
stant c s in Eq. ([TJ is independent of the closed-packed 
packing fraction 4> cp and the crystal structure C. 



3. Substitutionally disordered solid phase 

The density profiles of a crystal structure C with the 
sites occupied randomly by positively and negatively 
charged colloids with probability x\ and X2, respectively, 
are approximated, in place of Eq. ©, by 



i(r) =Xi^28(r-s). 



(8) 



sec 



Therefore, the Helmholtz free energy per particle is ap- 
proximated by 



c s + y Xi h\{xi) + \n(cf>) — 3 In 1 



(9) 



1/3- 



(2x-l) 2 ^ exp(-Ka(|s-s'|/a-2)) 



2K r T* 



(s,s')eC* 



s - s'\ a 



with C* = {(s,s') € C X C : s' € Wc,& ^ s'}. Note 
the "entropy of mixing" contribution due to the random 
occupation of sites and that the screened-Coulomb en- 
ergy vanishes at x = 1/2 for the present choice of equally 
sized, oppositely charged particles. 



4- Hard-sphere freezing 

In order to fix the values of the integration constants 
c f (see Eq. (4}) and c s (see Eqs. Q and ©) it is re- 
quired that hard-sphere freezing from a fluid to a ran- 
dom face-centered cubic (rfcc) phase, i.e., a substitution- 
ally disordered solid, takes place in the limit T* — >• oo 
at p* ocx = 11.54(4)tt/6 « 6.042 (see Ref. [13]). As the 
Gibbs free energies per particle of the fluid and the rfcc 
solid must be equal at coexistence, the two constants c/ 
and c s are related by Cf — c s « 1.084. As only differences 
of Gibbs free energies are relevant in order to determine 
the equilibrium state, we set c s — without restriction 
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TABLE I: Candidate solid structures chosen on the basis 
of Refs. 0, @-Q3]- For detailed structure information see 
Ref. [18fl . The second column indicates a crystal ("c") or 
a substitutionally disordered solid ("d"). The third and the 
fourth columns give the numbers of particles per unit cell K c 
and the closed-packed packing fractions (/> cp , respectively. The 
CuAu structure with the aspect ratio c/a G [1, y/2] of the un- 
derlying tetragonal lattice degenerates to the CsCl structure 
in the case c/a — 1. The LSg cc structure has been introduced 
in Ref. [ljj|. The NbP structure was called "tetragonal" in 
Ref. 0. 



of generality. The binodals of the fluid-rfcc solid coexis- 
tence region at high temperatures are found to be located 
at packing fractions 0f ~ 0.4916, r f cc ~ 0.5547, i.e., the 
coexistence region is in good agreement with that found 
in computer simulations [17j, albeit slightly wider. 



C. Phase diagrams 

In order to calculate phase diagrams in terms of T*, 
p* , and x a set of candidate solid structures is chosen in 
the next section. Let g(T* ,p* ,x) be defined as the mini- 
mum of the Gibbs free energies per particle constructed 
in the previous subsection of the fluid, the candidate 
crystals, and the candidate substitutionally disordered 
solids for given (T*,p*,x). The equilibrium Gibbs free 
energy is the convex hull of g(T* ,p* , x) as a function of 
x € [0, 1], which can be readily constructed numerically 
on a x-grid. During this calculation the coexistence re- 
gions of the phase diagrams can be identified as the set of 
state points where Maxwell's common tangent construc- 
tion applies. 

The Gibbs free energies per particle as defined in 
the previous subsection exhibit an unphysical feature at 
low pressures, where for all temperatures the fluid be- 
comes apparently unstable with respect to an rfcc crys- 
tal; the reason for this phenomenon is related to the 
well-known incorrect low-density asymptotics of the free- 
volume equation of state [lg|. In order to resolve this 
problem it is assumed that once a stable fluid state for 
given (T*,p*,x) is found, a stable fluid state is also as- 
sumed for (T* , a;) and all pressures smaller than p* . 
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FIG. 1: Phase diagrams for na = 3 at compositions x in terms of (T*jjj*) and {T*,cf>). The green dashed lines in panel (f) 
reproduce the phase diagram of the computer simulation study of Ref. [ljj]. Labels of pure solid phases correspond to Tab.|l] 
Panels (c) and (d) exhibit a two-phase coexistence region M2 = CsCl + rfcc. Horizontal thin lines in panels (b), (d), and (f) 
are tie lines. 



III. RESULTS AND DISCUSSION 

Table Q] lists the candidate solid structures considered 
here. This choice is based on the structures found in 



computer simulation studies of the cases x = (see 
Ref. [H) and na = 3,x = 1/2 (see Ref. [lOJ]), as well 
as on Refs. 0, El Bj, where the limit T* ->• is ad- 
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FIG. 2: Phase diagrams for na = 3 in terms of (p*,a;) at temperatures (a) T* = 10, (b) T* = 1, and (c) T* = 0.1 (see also 
Fig. [TJ. Horizontal thin lines in panels (b) and (c) are tie lines. 



dressed by means of Madelung energy sums. Moreover, 
the CsCl (cesium chloride), CuAu (copper gold), NaCl 
(sodium chloride), and C113AU structures have been iden- 
tified in experiments [3, If . Table U indicates whether 
the solid is crystalline ("c") or substitutionally disor- 
dered ( "d" ) and it exhibits the numbers of particles per 
unit cell K c as well as the closed-packed packing frac- 
tions 4> C p. Detailed structure information can be found 
in Ref. [la] . The CuAu structure, which is described by 
a tetragonal lattice of aspect ratio c/a e [1,^2] and a 
two-particle basis, degenerates to the CsCl structure in 
the case c/a = 1. The structure denoted by LS| CC was in- 
troduced in Ref. [l9j and the NbP (niobium phosphide) 
structure was called "tetragonal" in Ref. [10j . 



A. Case ko = 3 

The phase diagrams for the case Ka = 3 in terms of 
(T*,p*) and (T*,<f>) at compositions x = 0, x = 1/4, 
and x = 1/2 are displayed in Fig. [1] Thick solid lines 
represent phase boundaries whereas thin horizontal lines 
in (T*,cj)) diagrams (panels (b), (d), and (f)) are tie 
lines connecting coexisting states. At low temperatures 
two-phase coexistence Mi := CsCl + rfcc is found (see 
Figs. (He) and (d)). 

By construction a fluid-rfee transition takes place for 
any na and x in the limit T* — > 00 at the coexistence 
pressure p* oox and volume fractions <f>{ and 4> r f cc (see 
Sec. llTB~4|) . For T* -> and re = 1/2 a CsCl crystal 
coexists with a dilute gas, in agreement with Ref. [l2[ 
because the free energy of the present formalism reduces 
to Madelung-like energy sums in this limit. A cut of 
the phase diagrams at composition x — is shown in 
Figs. Ha) an d (b), which involves merely the fluid and 
the rfcc structures, which is consistent with computer 
simulation results [9| . Figures QJc) and (d) display a cut 
at x = 1/4, with an additional CU3AU phase as well 
as the two-phase coexistence region Mi- The structures 



at composition x = 1/4 differ from those at composi- 
tion x = 3/4 only by an exchange of the colloid species 
(lo2). 

In order to assess the reliability of the approach de- 
scribed in Sec. [TT1 Fig. QJf) compares the calculated 
phase diagram (red solid lines) in terms of (T* , 4>) for 
Ka = 3, x = 1/2 with that obtained by means of free en- 
ergy calculations using computer simulations in Ref. [10J] 
(green dashed lines) . Both studies agree in the predicted 
stable structures fluid, rfcc, CsCl, CuAu, and NbP. The 
overall topology is similar for both approaches; how- 
ever, the present formalism overestimates the stability 
of the NbP structure leading to a fiuid-NbP transition 
(see Figs.QJe) and (f)) which is not observed in the com- 
puter simulation. Moreover, our formalism leads to a 
fluid-CuAu transition in a very narrow window around 
T* fa 0.4 but no rfcc-CuAu transition is found, whereas 
it is the opposite situation with the computer simulation 
results. Agreement between mean- field theory and com- 
puter simulation is observed with respect to the order 
of the phase transitions: The CsCl-CuAu transition is of 
second order because the CuAu structure transforms con- 
tinuously into CsCl upon c/a — > 1, and the other phase 
transitions are of first order. Both the rfcc-NbP and the 
CuAu-NbP phase transitions are described as "weakly 
first-order" in Ref. [10], whereas within our approach, 
only the rfcc-NbP transition exhibits a very narrow but 
non-vanishing 0-gap and the CuAu-NbP transitions is 
strongly first-order (see Fig. QJf)). The quantitative dis- 
agreement in the strength of the first-order CuAu-NbP 
transition can be understood on the basis of a smearing 
out of structural differences due to fluctuations, which 
are present in computer simulations but which are not 
fully accounted for by mean-field theories. This com- 
parison between the formalism of Sec. [IT] and the com- 
puter simulation study of Ref. [1(| for the special case 
Ka = 3,x = 1/2 shows that, apart from well-known de- 
fects of mean-field theories, the present approach is semi- 
quantitatively reliable. Moreover, the simplicity of the 
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FIG. 3: Phase diagrams for na = 10 at compositions x in terms of (T* ,p*) and (T*, (j>). Labels of pure solid phases correspond 
to Tab.[U Panels (c) and (d) exhibit two-phase coexistence regions Mi — CsCl -(-fluid, Mi — CsCl + rfcc, and M$ — CuAu + rfcc. 
Horizontal thin lines in panels (b), (d), and (f) are tie lines. 



present formalism gives computational advantages over 
computer simulations such that now complete phase dia- 
grams in terms of (T* ,p*,x) as a function of the param- 
eter na can be determined readily. 



Figurc[5]displays the phase diagram for na — 3 in terms 
of (p*,x) for temperatures (a) T* = 10, (b) T* = 1, and 
(c) T* = 0.1. At high temperatures (see Fig. |2£a)) only 
fluid and rfcc structures are present and the fluid-rfcc 



7 



transition line becomes independent of the composition 
x in the limit T* — > oo. At lower temperatures (see 
Figs. GOT.) and (c)) CsCl, CuAu, Cu 3 Au, and NbP crystal 
structures occur at fixed compositions. At low tempera- 
tures and pressures as well as strongly asymmetric mix- 
tures (see Fig. [2jc) for p* < 1 and x < 1/4 ori> 3/4) 
the MSA approximation applied to model the fluid phase 
leads to an unphysical artifact which exhibits the appar- 
ent coexistence of an almost pure rfcc crystal with a less 
pure fluid. The reason for this unphysical phenomenon 
is that under these conditions the MSA pair distribu- 
tion function becomes negative such that an increasing 
repulsive interaction potential leads to a more and more 
negative, i.e., attractive, contribution to the free energy. 
However, outside of this range of the phase diagram MSA 
leads to a physically reasonable description of the fluid 
phase. The full phase diagram for na = 3 in terms of 
(T* , p* , x) can be inferred from the two-dimensional cuts 
in Figs. Q] and H 



B. 



Case na = 10 



In order to study the changes of the phase diagram 
upon changing the screening strength na, the phase di- 
agrams for various screening strengths na £ {2, 3, 5, 10} 
have been calculated. The trends in the variations of the 
phase diagrams upon changing na are found to be mono- 
tonic such that it is sufficient in the following to consider 
the case na = 10. 

Figures [3] displays the phase diagram for na — 10 in 
terms of (T*,p*) and (T*,<f>) at compositions x = 0, 
x = 1/4, and x — 1/2. Two-phase coexistence re- 
gions M x \= CsCl + fluid, M 2 = CsCl + rfcc, and 
M3 := CuAu + rfcc are present in Figs.^c) and (d). By 
comparison of Figs. Q] and [3] one infers a shift of the rfcc 
phase to lower temperatures upon increasing na. This 
observation can be understood by the fact that the inter- 
action potential is approaching the hard-sphere potential 
in the limit na — ¥ 00. The CU3AU phase in Figs. [TJc) 
and [3fc) or Figs. QJd) and [3Jd) shrinks upon increas- 
ing na, which is partly due to the growing rfcc phase. 
Figures QJe) an d EJe) or Figs. QJf ) and [SJf) exhibit an 
increasing temperature range of stability of the CuAu 
structure upon increasing na. As a consequence the NbP 
phase, located in between the extending rfcc and CuAu 
phases, shrinks upon increasing na. Moreover, upon in- 
creasing na, the fluid-NbP transition disappears and an 
rfcc-CuAu coexistence is established. Finally, the <j> range 
of the CsCl phase becomes smaller upon increasing na as 
can be inferred from Figs.QTf) and^f). 

In order to make predictions on the conditions to 
synthesize certain crystal structures, it is expected 
within the present formalism that colloids with strongly 



screened Coulomb interaction, i.e., large values of na, 
are preferable to prepare CuAu structures, whereas CsCl, 
CU3AU, and NbP crystals are expected to be found most 
easily in systems of weakly screened Coulomb interaction. 
Given a certain crystal structure has been prepared, the 
above reasoning leads to the following conclusions: CuAu 
structures become more whereas CsCl, CU3AU, and NbP 
structures become less stable against temperature varia- 
tions upon increasing na. 



IV. SUMMARY 

In this work a simple, yet semi-quantitative mean-field 
theory for binary colloidal mixtures of spherical particles 
of equal radii and opposite charges has been described 
and applied to calculate the global phase diagrams in 
terms of temperature, pressure and composition for var- 
ious screening strengths (see Figs. Q~r|3]). The formalism 
interpolates between the hard-sphere limit at high tem- 
peratures and a Madehmg-like description at low temper- 
atures. The reliability of the method has been checked 
by comparison with computer simulation results [9|, [l0( 
for one single composition and screening strength as well 
as with Madelung sum calculations 0, LLJ], [l2[ in the 
low-temperature limit. Limitations of the present for- 
malism can be understood as a result of the absence of 
fluctuations within mean-field theories and properties of 
the mean-spherical approximation (MSA). The influence 
of the screening strength on the stability of crystalline 
phases has been discussed, which has implications, e.g., 
on the preparation procedure and the choice of experi- 
mental conditions. 

An extension of the presented theory to multi- 
component, size- and charge-polydisperse colloidal sus- 
pensions along the lines described in Sec. [TT] is straight- 
forward as the free volume equation of state used for 
the solid phases does not explicitly depend on the col- 
loid radii and more general formulations of the MSA are 
available [2Cj, [21| . It is therefore a matter of choosing an 
appropriate set of candidate solid structures (see Tab. [I]), 
which could be motivated by findings from experiments 
and computer simulations of special cases. 
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